home *** CD-ROM | disk | FTP | other *** search
/ Graphics Plus / Graphics Plus.iso / general / procssng / ccs / ccs-11tl.lha / lbl / hips / sources / tools / mhisto.c < prev    next >
Encoding:
C/C++ Source or Header  |  1992-06-06  |  7.7 KB  |  295 lines
  1. /*
  2. % mhisto.c - compute grey-level histograms for all formatted images.
  3. %
  4. %    Copyright (c)    1991    Jin, Guojun
  5. %
  6. % usage:
  7. %    mhisto [-l(#logbins)] [-c] [-m #] [-b # -n #] [-r] [-f -t -z] [-s #]
  8. %        [-a #] [-S [#]] [-v [#]] [-M] [<] image [> [-o] outhist]
  9. % where log(bins) is the log2 of the number of bins. Default bins is 256.    */
  10. char    usage[]="options\n\
  11. -c causes multiple frame sequences to collapse to a single histogram,\n\
  12.    instead of a separate histogram being generated for each input frame.\n\
  13.     It can be used to compute 3D image histogram.\n\
  14. -l #    bin # in log2 (default=8)\n\
  15. -m #    specify the max output value. The maximum integer input is 2^31.\n\
  16. -1 -2 -a #    any file (1=short, 2=int) with size #\n\
  17. -b #    begin process from #th frame.\n\
  18. -n #    to process # frames.\n\
  19. -f -t    will eliminate the frequent top or zero value count.    \n\
  20. -r    recalculate maximum count for each frame and -R option will be used\n\
  21.     in mdisphist to Retrieve these conuts.\n\
  22. -s #    set display scale (maximum count) for entire histogram.    \n\
  23. -v [#]    dig valley at given position [#], which defaulted at 32, and the valley\n\
  24.     is digged to 2/3 at position left and to 1/3 at its right.\n\
  25. -V #    Valley width. Default = 9.\n\
  26. -z    count zeros. Default is eliminate zero value.    \n\
  27. -M    allow to display some important message.    \n\
  28. -S [#]    smooth the histogram. The smaller following number, the rougher of\n\
  29.     smoothing.    \n\
  30. [<] image [> histogram]    \n";
  31. /*
  32. % compile:    cc -O -o mhisto mhisto.c -lscs3 -lccs -lhips -lrle -ltiff -lm
  33. % note:    there are 2 integers (bin_width and bins) between header and data and
  34. %    1 integer (maxcnt) before each frame.
  35. %
  36. % AUTHOR:    Jin Guojun - LBL    1/8/91
  37. */
  38.  
  39. #include "header.def"
  40. #include "imagedef.h"
  41. #include <math.h>
  42.  
  43. bool    any_pure, Msg, cflag, zflag, recalc, topf, freqf, set_f, smooth;
  44. U_IMAGE    uimg;
  45.  
  46. #define    inbuf    uimg.src
  47. #define    frame    uimg.frames
  48.  
  49. #ifndef    MaxBin
  50. #define    MaxBin  131072
  51. #endif
  52. #ifndef    SMOOTH
  53. #define    SMOOTH    17
  54. #endif
  55.  
  56. struct    {
  57.     MType    max, frm, bin;
  58.     }rec;
  59.  
  60. struct    {
  61.     int    pos, w;
  62.     }valley={0, 9};
  63.  
  64. #define    GValue()    arget(ac, av, &f, &fr)
  65.  
  66.  
  67. main(ac, av)
  68. int    ac;
  69. char**    av;
  70. {
  71. MType    f, fr, nfrms=0, bgnf=0, MaxVal=256,
  72.     maxcnt=0, secmax, maxpos, secmaxp, numbin=0,
  73.     binwidth,/* reduce resolution. i.e. if binwidth=3, then 1,2,3
  74.             will go to same bin and 4,5,6 will go to same bin */
  75.     fsize, simu_size=0;
  76. register MType    *hist;
  77. float    scale;
  78.  
  79. format_init(&uimg, IMAGE_INIT_TYPE, HIPS, -1, *av, "S20-1");
  80.  
  81. for (f=1; f<ac; f++)
  82.     if (*av[f] == '-')    {
  83.     fr = 1;    
  84.     switch(av[f][fr++])    {
  85.     case 'D':debug++;    break;
  86.     case 'M':Msg++;    break;
  87.     case 'S':
  88.         smooth = GValue();
  89.         if (smooth<1 || smooth>SMOOTH)    smooth=5;
  90.         break;
  91.     case '1': case '2':    simu_size = av[f][fr-1] - '0';
  92.     case 'a':    any_pure = GValue();    break;
  93.     case 'b':
  94.         bgnf = GValue() - 1;    break;
  95.     case 'c':
  96.         cflag++;    break;
  97.     case 'f':    freqf++;    break;
  98.     case 'l':
  99.         numbin = 1 << (int)GValue();    break;
  100.     case 'm':
  101.         MaxVal = GValue();    break;
  102.     case 'n':
  103.         nfrms = GValue();    break;
  104.     case 'r':
  105.         recalc++;    break;
  106.     case 's':
  107.         set_f = maxcnt = GValue();    break;
  108.     case 't':
  109.         topf++;    break;
  110.     case 'v':
  111.         valley.pos = GValue();
  112.         if (!valley.pos)    valley.pos = 32;
  113.         break;
  114.     case 'V':
  115.         valley.w = GValue();    break;
  116.     case 'z':
  117.         zflag++;    break;
  118.     case 'o':
  119.         if (out_fp=fopen(av[f]+2, "wb"))    break;
  120.     default:
  121. info:        usage_n_options(usage, f, av[f]);
  122.     }
  123.     }
  124.     else if (!(in_fp=fopen(uimg.name=av[f], "rb")))
  125.         syserr("%s - not found", av[f]);
  126.  
  127. io_test(fileno(in_fp), goto    info);
  128.  
  129. if (any_pure > 0)    {
  130.     fsize = uimg.width = any_pure>>simu_size;
  131.     f = frame = uimg.height = 1;
  132.     uimg.pxl_in = 1 << (uimg.in_form=simu_size);
  133.     (*uimg.std_swif)(FI_INIT_NAME, &uimg, uimg.name, 0);
  134. }
  135. else    {
  136.     (*uimg.header_handle)(HEADER_READ, &uimg, 0, 0);
  137.     fsize = uimg.width * uimg.height;
  138.     f = frame;
  139.     if (bgnf<0)    bgnf=0;
  140.     else if(bgnf>=frame)    bgnf = frame-1;
  141.     if (bgnf)    fr = fseek(in_fp, fsize*bgnf, 1);
  142.     if (fr==EOF)    syserr("passing %d first %d frames", fr, bgnf);
  143.     if (Msg)    msg("fp at %d\n", fr);
  144.     if (nfrms<1 || nfrms > frame-bgnf)    nfrms=frame-bgnf;
  145.     if (nfrms)    f = frame = nfrms;/*    display less frames    */
  146. }
  147. if (uimg.in_form != IFMT_BYTE)
  148.     MaxVal = 65536;
  149. if (!numbin)
  150.     numbin = 256;
  151. else if (numbin > MaxVal)
  152.     numbin = MaxVal;
  153.  
  154. if (uimg.in_form > IFMT_FLOAT)
  155.     syserr("image format must be lower than float point");
  156.  
  157. hist = zalloc(numbin, (MType)sizeof(*hist), "hist");
  158.  
  159. uimg.o_form = IFMT_HIST;
  160. uimg.pxl_out = sizeof(int);
  161.  
  162. if (cflag)    frame = 1;
  163.  
  164. binwidth = MaxVal / numbin;
  165. message("maxv=%ld, numbin=%ld, bin_w=%ld, fs=%ld",MaxVal,numbin,binwidth,fsize);
  166. if (valley.pos)
  167.     message(", valley->%d(width=%d)", valley.pos-(valley.w<<1)/3, valley.w);
  168. mesg("\n");
  169.  
  170. (*uimg.header_handle)(HEADER_WRITE, &uimg, ac, av, True);
  171.  
  172. fwrite(&binwidth, 1, sizeof(binwidth), out_fp);
  173. fwrite(&numbin, 1, sizeof(numbin), out_fp);
  174.  
  175. for(fr=0; fr<f; fr++) {
  176. register MType    i;
  177.     (*uimg.std_swif)(FI_LOAD_FILE, &uimg, uimg.load_all=0, No);
  178.     if (recalc)    maxcnt = 1;
  179.     switch(uimg.in_form)
  180.     {
  181.     case IFMT_BYTE:{
  182.         register byte* bp=inbuf;
  183.         register MType    j,k;
  184.             for (i=0; i<fsize; i++) {
  185.                 j = *bp++ & 0xFF;
  186.                 if (!(j | zflag))    continue;
  187.                 k = ++hist[j/binwidth];
  188.                 if (k>maxcnt && k>set_f){
  189.                    secmaxp = maxpos;    maxpos = j;
  190.                    if (set_f){
  191.                      secmax = set_f;    set_f = k;
  192.                    }
  193.                    else{ secmax = maxcnt;    maxcnt = k; }
  194.                 }
  195.             }
  196.         }
  197.         break;
  198.     case IFMT_LONG: {
  199.         register int*    ip=inbuf, min=16777276, max=0;
  200.         register unsigned short    *bp=inbuf;
  201.             for (i=0; i<fsize; i++, ip++)
  202.                 if (*ip < min)    min = *ip;
  203.                 else if (*ip > max)    max = *ip;
  204.             i = abs(max - min);
  205.             scale = (float)MaxVal / i;
  206.     if(Msg)    message("min=%d, max=%d, scale=%f\n", min, max, scale);
  207.             ip = inbuf;
  208.             for (i=0; i<fsize; i++)
  209.                 *bp++ = (*ip++ - min) * scale + .5;
  210.         }
  211.         goto    Hshort;
  212.     case IFMT_FLOAT:{
  213.         register float*    fp=inbuf, min=1e78, max = -1e78;
  214.         register unsigned short    *bp=inbuf;
  215.             for (i=0; i<fsize; i++, fp++)
  216.                 if (*fp < min)    min = *fp;
  217.                 else if (*fp > max)    max = *fp;
  218.             i = max - min;
  219.             scale = (float)MaxVal / i;
  220.     if(Msg)    message("min=%f, max=%f, scale=%f\n", min, max, scale);
  221.             fp = inbuf;
  222.             for (i=0; i<fsize; i++)
  223.                 *bp++ = (*fp++ - min) * scale + .5;
  224.         }
  225.     case IFMT_SHORT:
  226. Hshort:        {
  227.         register unsigned short    *sp=inbuf;
  228.         register MType    j,k;
  229.             for (i=0; i<fsize; i++) {
  230.                 j = *sp++;
  231.                 if (!(j | zflag))    continue;
  232.                 k = ++hist[j/binwidth];
  233.                 if (k>maxcnt && k>set_f){
  234.                    secmaxp = maxpos;    maxpos = j;
  235.                    if (set_f){
  236.                      secmax = set_f;    set_f = k;
  237.                    }
  238.                    else{ secmax = maxcnt;    maxcnt = k; }
  239.                 }
  240.             }
  241.         }
  242.         break;
  243.     }
  244.     if (smooth){
  245.     register int    j, s, n = (smooth << 1) + 1;
  246.         for (i=smooth; i<numbin-smooth; i++){
  247.             for (s=0, j = -smooth; j<=smooth; j++)
  248.                 s += hist[i+j];
  249.             hist[i] = s / n;
  250.         }
  251.     }
  252.     if (valley.pos){
  253.     register int    i, p=valley.pos - (valley.w<<1)/3;
  254.         for (i=valley.w; i; i--)
  255.             hist[p++] = 0;
  256.     }
  257.     if (topf)    hist[numbin-1] = 0;
  258.     if (freqf && maxpos && maxpos != numbin-1){
  259.         hist[maxpos]=0;
  260.         if (set_f)    set_f=secmax;
  261.         else    maxcnt = secmax;
  262.         maxpos = secmaxp;
  263.     }
  264.     if (set_f > rec.max){
  265.         rec.max = set_f;
  266.         rec.frm = fr;    rec.bin = maxpos;
  267.     }
  268.     else if (!set_f && maxcnt > rec.max){
  269.         rec.max = maxcnt;
  270.         rec.frm = fr;
  271.         rec.bin = maxpos;
  272.     }
  273.     if(Msg)    {
  274.     register int    maxc = set_f ? set_f : maxcnt;
  275.         message("%s:(%03d) max count = %d at postion %d\n",
  276.             *av, fr+1, maxc, maxpos);
  277.     }
  278.     if (!cflag) {
  279.         if (fwrite(&maxcnt, 1, sizeof(maxcnt), out_fp) != sizeof(maxcnt))
  280.             syserr("can not write max count");
  281.         if (fwrite(hist, sizeof(numbin), numbin, out_fp) != numbin)
  282.                 syserr("error during write");
  283.         for (i=0; i<numbin; i++)    hist[i] = 0;
  284.     }
  285.     }
  286. if (cflag){
  287.     if (fwrite(&maxcnt, 1, sizeof(maxcnt), out_fp) != sizeof(maxcnt))
  288.         syserr("can not write max count");
  289.     if (fwrite(hist, sizeof(*hist), numbin, out_fp) != numbin)
  290.         syserr("error during write");
  291. }
  292. message("%s: max count = %ld at frame %ld, bin %ld\n",
  293.     *av, rec.max, rec.frm, rec.bin);
  294. }
  295.